Human CTB BPS Project
Sean Nguyen
Overall Objective
Load Libraries
library(tidyverse)
library(cowplot)
library(broom)
library(plotly)Import data
Convert data from ‘wide’ to ‘long’ format
# data
data1 <- data %>%
gather(Sample,Count,2:49)
# Separate samples by identifiers
data2 <- data1 %>%
separate(Sample, into=c("Sample_ID","Dilution_factor",
"Injection","Tech_rep", sep = "_")) %>%
select(-`_`)
# Standards
standards1 <- standards %>%
gather(Sample,Count,2:25)
standards2 <- standards1 %>%
separate(Sample, into=c("Sample_ID","When","Dilution_factor",
"Nano_day","Injection","Tech_Rep", sep = "_")) %>%
select(-`_`)Factor the data into categorical variables
# Refactoring Columns for samples
data2$Sample_ID <- as.factor(data2$Sample_ID)
data2$Dilution_factor <- as.numeric(data2$Dilution_factor)
data2$Injection<- as.factor(data2$Injection)
data2$Tech_rep <- as.numeric(data2$Tech_rep)
data2## # A tibble: 48,000 × 6
## particle_size Sample_ID Dilution_factor Injection Tech_rep Count
## * <dbl> <fctr> <dbl> <fctr> <dbl> <int>
## 1 0.5 1 25 1 0 0
## 2 1.5 1 25 1 0 0
## 3 2.5 1 25 1 0 0
## 4 3.5 1 25 1 0 0
## 5 4.5 1 25 1 0 0
## 6 5.5 1 25 1 0 0
## 7 6.5 1 25 1 0 0
## 8 7.5 1 25 1 0 0
## 9 8.5 1 25 1 0 0
## 10 9.5 1 25 1 0 0
## # ... with 47,990 more rows# Refactoring COlumns for key
key$Sample_ID <- as.factor(key$Sample_ID)
key$Time <- as.factor(key$Time)
key$Treatment <- as.factor(key$Treatment)
key$Volume <- as.numeric(key$Volume)
key$Treatment <- factor(key$Treatment,levels = c('DMSO','EGF','BPS','BPS_EGF'))
key## # A tibble: 8 × 4
## Sample_ID Time Treatment Volume
## <fctr> <fctr> <fctr> <dbl>
## 1 1 48 BPS 145
## 2 2 48 DMSO 112
## 3 3 48 BPS_EGF 110
## 4 4 48 EGF 120
## 5 5 96 DMSO 105
## 6 6 96 BPS_EGF 112
## 7 7 96 EGF 105
## 8 8 96 BPS 110# Refactoring columns for standards
standards2$Sample_ID <- as.factor(standards2$Sample_ID)
standards2$When <- as.factor(standards2$When)
standards2$Dilution_factor <- as.numeric(standards2$Dilution_factor)
standards2$Injection <- as.factor(standards2$Injection)
standards2$Nano_day <- as.numeric(standards2$Nano_day)
standards2## # A tibble: 24,000 × 8
## particle_size Sample_ID When Dilution_factor Nano_day Injection
## * <dbl> <fctr> <fctr> <dbl> <dbl> <fctr>
## 1 0.5 std after 125 1 1
## 2 1.5 std after 125 1 1
## 3 2.5 std after 125 1 1
## 4 3.5 std after 125 1 1
## 5 4.5 std after 125 1 1
## 6 5.5 std after 125 1 1
## 7 6.5 std after 125 1 1
## 8 7.5 std after 125 1 1
## 9 8.5 std after 125 1 1
## 10 9.5 std after 125 1 1
## # ... with 23,990 more rows, and 2 more variables: Tech_Rep <chr>,
## # Count <int>Back calculate standards
standards2 <- standards2 %>%
mutate(True_Count=Dilution_factor*Count)
# Set the correct order of 'categorical factors'
standards2$Nano_day <- factor(standards2$Nano_day, levels=c('1','2'))
standards2$When <- factor(standards2$When, levels=c('before','after'))
standards2## # A tibble: 24,000 × 9
## particle_size Sample_ID When Dilution_factor Nano_day Injection
## <dbl> <fctr> <fctr> <dbl> <fctr> <fctr>
## 1 0.5 std after 125 1 1
## 2 1.5 std after 125 1 1
## 3 2.5 std after 125 1 1
## 4 3.5 std after 125 1 1
## 5 4.5 std after 125 1 1
## 6 5.5 std after 125 1 1
## 7 6.5 std after 125 1 1
## 8 7.5 std after 125 1 1
## 9 8.5 std after 125 1 1
## 10 9.5 std after 125 1 1
## # ... with 23,990 more rows, and 3 more variables: Tech_Rep <chr>,
## # Count <int>, True_Count <dbl>Summarize three technical standard replicates
standards3 <- standards2 %>%
group_by(particle_size,Sample_ID,When,Dilution_factor,Nano_day,Injection) %>%
summarise( tech_N = length(True_Count),
tech_mean = mean(True_Count),
tech_sd = sd(True_Count),
tech_se = tech_sd/sqrt(tech_N))
standards3## Source: local data frame [8,000 x 10]
## Groups: particle_size, Sample_ID, When, Dilution_factor, Nano_day [?]
##
## particle_size Sample_ID When Dilution_factor Nano_day Injection
## <dbl> <fctr> <fctr> <dbl> <fctr> <fctr>
## 1 0.5 std before 125 1 1
## 2 0.5 std before 125 1 2
## 3 0.5 std before 125 2 1
## 4 0.5 std before 125 2 2
## 5 0.5 std after 125 1 1
## 6 0.5 std after 125 1 2
## 7 0.5 std after 125 2 1
## 8 0.5 std after 125 2 2
## 9 1.5 std before 125 1 1
## 10 1.5 std before 125 1 2
## # ... with 7,990 more rows, and 4 more variables: tech_N <int>,
## # tech_mean <dbl>, tech_sd <dbl>, tech_se <dbl>Summarize standards by injection
standards4 <- standards3 %>%
group_by(Nano_day,When,particle_size) %>%
summarise( inj_N = length(tech_mean),
inj_mean = mean(tech_mean),
inj_sd = sd(tech_mean),
inj_se = inj_sd/sqrt(inj_N))
standards4## Source: local data frame [4,000 x 7]
## Groups: Nano_day, When [?]
##
## Nano_day When particle_size inj_N inj_mean inj_sd inj_se
## <fctr> <fctr> <dbl> <int> <dbl> <dbl> <dbl>
## 1 1 before 0.5 2 0 0 0
## 2 1 before 1.5 2 0 0 0
## 3 1 before 2.5 2 0 0 0
## 4 1 before 3.5 2 0 0 0
## 5 1 before 4.5 2 0 0 0
## 6 1 before 5.5 2 0 0 0
## 7 1 before 6.5 2 0 0 0
## 8 1 before 7.5 2 0 0 0
## 9 1 before 8.5 2 0 0 0
## 10 1 before 9.5 2 0 0 0
## # ... with 3,990 more rowsPlot before and after plots, facet by experimental day
std_plot <- standards4 %>%
ggplot(aes(x = particle_size, y = inj_mean, color=When))+
geom_line(size=2) + xlim(0,300)+ #line size, x-axis scale
geom_ribbon(aes(ymin=inj_mean-inj_se, ymax=inj_mean+inj_se),
alpha=0.4,fill = alpha('grey12', 0.2)) + #error bars
scale_y_continuous(expand=c(0,0))+ #set bottom of graph
xlab("Particle Size") + # X axis label
ylab("\nMean Particle Concentration/ml\n") + # Y axis label
ggtitle("Nanosight Histogram of\n100nm Standards")+ #title
labs(color="Condition")+ #Label table title
facet_wrap(~ Nano_day)
std_plot## Warning: Removed 1400 rows containing missing values (geom_path).
# ggsave("Standards_histogram_plot.png",
# height = 5, width = 7, dpi = 300, units= "in")Standards particle concentrations from each experimental day
standards_df <- standards4 %>%
group_by(Nano_day,When) %>%
summarise(total=sum(inj_mean))
standards_df## Source: local data frame [4 x 3]
## Groups: Nano_day [?]
##
## Nano_day When total
## <fctr> <fctr> <dbl>
## 1 1 before 45989600958
## 2 1 after 41407513229
## 3 2 before 51412158250
## 4 2 after 50587138812Bar graph of standards particle concentrations
standards_bar <- standards_df %>%
ggplot(aes(x=Nano_day,y=total,fill=When))+
geom_col(position="dodge")+
scale_y_continuous(expand=c(0,0))+ #set bottom of graph
xlab("Experimental Day") + # X axis label
ylab("\nMean Particle Concentration/ml\n") + # Y axis label
ggtitle("Nanosight Histogram of\n100nm Standards")+ #title
labs(color="When") #Label table title
standards_bar
# ggsave("Standards_bar_plot.png",
# height = 5, width = 7, dpi = 300, units= "in")Intraassay variability
Intra.assay_cv <- standards_df %>%
group_by(Nano_day) %>%
summarise(Intra_Day_N = length(total),
Intra_Day_mean = mean(total),
Intra_Day_sd = sd(total),
Intra_Day_se = Intra_Day_sd/sqrt(Intra_Day_N),
Intra_Day_cv = Intra_Day_sd/Intra_Day_mean )
Intra.assay_cv## # A tibble: 2 × 6
## Nano_day Intra_Day_N Intra_Day_mean Intra_Day_sd Intra_Day_se
## <fctr> <int> <dbl> <dbl> <dbl>
## 1 1 2 43698557094 3240025305 2291043865
## 2 2 2 50999648531 583376839 412509719
## # ... with 1 more variables: Intra_Day_cv <dbl># # Save as .csv
# write_csv(Intra.assay_cv,"Intra.assay_cv.csv")Interassay variability
Inter.assay_cv <- Intra.assay_cv %>%
summarise(Inter_Day_N = length(Intra_Day_mean),
Inter_Day_mean = mean(Intra_Day_mean),
Inter_Day_sd = sd(Intra_Day_mean),
Inter_Day_se = Inter_Day_sd/sqrt(Inter_Day_N),
Inter_Day_cv = Inter_Day_sd/Inter_Day_mean )
Inter.assay_cv## # A tibble: 1 × 5
## Inter_Day_N Inter_Day_mean Inter_Day_sd Inter_Day_se Inter_Day_cv
## <int> <dbl> <dbl> <dbl> <dbl>
## 1 2 47349102812 5162651266 3650545719 0.1090338# # Save as .csv
# write_csv(Inter.assay_cv,"Inter.assay_cv.csv")Sample analysis
Back calculate the original concentration of the sample
data2 <- data2 %>%
mutate(True_Count = Dilution_factor*Count)
data2## # A tibble: 48,000 × 7
## particle_size Sample_ID Dilution_factor Injection Tech_rep Count
## <dbl> <fctr> <dbl> <fctr> <dbl> <int>
## 1 0.5 1 25 1 0 0
## 2 1.5 1 25 1 0 0
## 3 2.5 1 25 1 0 0
## 4 3.5 1 25 1 0 0
## 5 4.5 1 25 1 0 0
## 6 5.5 1 25 1 0 0
## 7 6.5 1 25 1 0 0
## 8 7.5 1 25 1 0 0
## 9 8.5 1 25 1 0 0
## 10 9.5 1 25 1 0 0
## # ... with 47,990 more rows, and 1 more variables: True_Count <dbl>Average three technical readings
data3 <- data2 %>%
group_by(particle_size,Sample_ID,Dilution_factor,Injection) %>%
summarise( tech_N = length(True_Count),
tech_mean = mean(True_Count),
tech_sd = sd(True_Count),
tech_se = tech_sd/sqrt(tech_N))
data3## Source: local data frame [16,000 x 8]
## Groups: particle_size, Sample_ID, Dilution_factor [?]
##
## particle_size Sample_ID Dilution_factor Injection tech_N tech_mean
## <dbl> <fctr> <dbl> <fctr> <int> <dbl>
## 1 0.5 1 25 1 3 0
## 2 0.5 1 25 2 3 0
## 3 0.5 2 25 1 3 0
## 4 0.5 2 25 2 3 0
## 5 0.5 3 25 1 3 0
## 6 0.5 3 25 2 3 0
## 7 0.5 4 25 1 3 0
## 8 0.5 4 25 2 3 0
## 9 0.5 5 25 1 3 0
## 10 0.5 5 25 2 3 0
## # ... with 15,990 more rows, and 2 more variables: tech_sd <dbl>,
## # tech_se <dbl>Summarize samples by injection (average both injections)
data4 <- data3 %>%
group_by(particle_size,Sample_ID,Dilution_factor) %>%
summarise( inj_N = length(tech_mean),
inj_mean = mean(tech_mean),
inj_sd = sd(tech_mean),
inj_se = inj_sd/sqrt(inj_N))
data4## Source: local data frame [8,000 x 7]
## Groups: particle_size, Sample_ID [?]
##
## particle_size Sample_ID Dilution_factor inj_N inj_mean inj_sd inj_se
## <dbl> <fctr> <dbl> <int> <dbl> <dbl> <dbl>
## 1 0.5 1 25 2 0 0 0
## 2 0.5 2 25 2 0 0 0
## 3 0.5 3 25 2 0 0 0
## 4 0.5 4 25 2 0 0 0
## 5 0.5 5 25 2 0 0 0
## 6 0.5 6 25 2 0 0 0
## 7 0.5 7 25 2 0 0 0
## 8 0.5 8 25 2 0 0 0
## 9 1.5 1 25 2 0 0 0
## 10 1.5 2 25 2 0 0 0
## # ... with 7,990 more rows# Average technical replicates and merge with key
merge <- left_join(key,data3, by= "Sample_ID")
merge## # A tibble: 16,000 × 11
## Sample_ID Time Treatment Volume particle_size Dilution_factor
## <fctr> <fctr> <fctr> <dbl> <dbl> <dbl>
## 1 1 48 BPS 145 0.5 25
## 2 1 48 BPS 145 0.5 25
## 3 1 48 BPS 145 1.5 25
## 4 1 48 BPS 145 1.5 25
## 5 1 48 BPS 145 2.5 25
## 6 1 48 BPS 145 2.5 25
## 7 1 48 BPS 145 3.5 25
## 8 1 48 BPS 145 3.5 25
## 9 1 48 BPS 145 4.5 25
## 10 1 48 BPS 145 4.5 25
## # ... with 15,990 more rows, and 5 more variables: Injection <fctr>,
## # tech_N <int>, tech_mean <dbl>, tech_sd <dbl>, tech_se <dbl># Save as .csv
# write_csv(merge,"Technical_replicate_average.csv")
# Average injection replicates and merge with key
merge1 <- left_join(key,data4, by= "Sample_ID")
merge1## # A tibble: 8,000 × 10
## Sample_ID Time Treatment Volume particle_size Dilution_factor inj_N
## <fctr> <fctr> <fctr> <dbl> <dbl> <dbl> <int>
## 1 1 48 BPS 145 0.5 25 2
## 2 1 48 BPS 145 1.5 25 2
## 3 1 48 BPS 145 2.5 25 2
## 4 1 48 BPS 145 3.5 25 2
## 5 1 48 BPS 145 4.5 25 2
## 6 1 48 BPS 145 5.5 25 2
## 7 1 48 BPS 145 6.5 25 2
## 8 1 48 BPS 145 7.5 25 2
## 9 1 48 BPS 145 8.5 25 2
## 10 1 48 BPS 145 9.5 25 2
## # ... with 7,990 more rows, and 3 more variables: inj_mean <dbl>,
## # inj_sd <dbl>, inj_se <dbl># #Save as .csv
# write_csv(merge1,"Injection_replicate_average.csv")Quick visualizations
Graphing all samples
sample_plot <- merge %>%
ggplot(aes(x=particle_size, y=tech_mean,color=Injection ))+ #plot
geom_ribbon(aes(ymin=tech_mean-tech_se,
ymax=tech_mean+tech_se),
alpha=0.2,fill = alpha('grey12', 0.2)) + #error bars
geom_line(size=2.0, alpha = 0.8) + xlim(0,500)+ #line size, x-axis scale
scale_y_continuous(expand=c(0,0))+ #set bottom of graph
xlab("Particle Size") + # X axis label
ylab("\nMean Particle Concentration/ml\n") + # Y axis label
ggtitle("Nanosight Histogram of\nhCTBS treated with BPS")+ #title
labs(color="Injection")+ #Label table title
facet_grid(Time ~ Treatment)
# geom_vline(xintercept = 200)+
# annotate("text", x= 350, y = 1E8, label= "200nm")
sample_plot
# ggsave("Nanosight_Sample_Histogram.png", plot = sample_plot,
# height = 10, width = 14, dpi = 200, units= "in")Interactive Plot
ggplotly(sample_plot)Particle concentration values for each of the samples
merge2 <- merge1 %>%
group_by(Time, Treatment, Volume) %>%
summarise(particle_conc=sum(inj_mean))
merge2## Source: local data frame [8 x 4]
## Groups: Time, Treatment [?]
##
## Time Treatment Volume particle_conc
## <fctr> <fctr> <dbl> <dbl>
## 1 48 DMSO 112 5332908929
## 2 48 EGF 120 13941558392
## 3 48 BPS 145 5733223262
## 4 48 BPS_EGF 110 5658288538
## 5 96 DMSO 105 5570460917
## 6 96 EGF 105 4416285746
## 7 96 BPS 110 5557785942
## 8 96 BPS_EGF 112 4699700962Correct for resuspension volume
merge3 <- merge2 %>%
mutate(particle_count = (Volume/1000)*particle_conc, # Create new column with number of particles
corrected_particle_conc = (particle_conc/1E9)) # Create new column with correct particle concentration
merge3## Source: local data frame [8 x 6]
## Groups: Time, Treatment [8]
##
## Time Treatment Volume particle_conc particle_count
## <fctr> <fctr> <dbl> <dbl> <dbl>
## 1 48 DMSO 112 5332908929 597285800
## 2 48 EGF 120 13941558392 1672987007
## 3 48 BPS 145 5733223262 831317373
## 4 48 BPS_EGF 110 5658288538 622411739
## 5 96 DMSO 105 5570460917 584898396
## 6 96 EGF 105 4416285746 463710003
## 7 96 BPS 110 5557785942 611356454
## 8 96 BPS_EGF 112 4699700962 526366508
## # ... with 1 more variables: corrected_particle_conc <dbl># Save as .csv
# write_csv(merge3,"Adjusted_particle_concentration.csv")Barplot
plot1 <- merge3 %>%
ggplot(aes(x = Treatment, y = corrected_particle_conc, fill = Treatment)) +
geom_bar(aes(text = paste("Particle Concentration:",
corrected_particle_conc)),
stat="identity", position = "dodge")+
xlab("\nTreatment\n") + # X axis label
ylab("\nMean Vessicle Concentration\n(10^9 particles/ ml)\n") + # Y axis label
ggtitle("Effect of BPS on Extracellular Vessicle\nRelease of hCTBs\n")+
scale_y_continuous(breaks = seq(0,14,2),
limits = c(0,14),
expand = c(0,0))+ # set bottom of graph
labs(color="Condition")+ # Label table title
facet_wrap(~Time)## Warning: Ignoring unknown aesthetics: textplot1
# ggsave("BPS_treated_hCTBs_48_96_facet_plot.png",
# height = 8, width = 10, dpi = 600, units= "in")Interactive Plot
ggplotly(plot1)Barplot
plot2 <- merge3 %>%
group_by(Time) %>%
ggplot(aes(x = Treatment, y = corrected_particle_conc, fill = Time )) +
geom_bar(position = "dodge", stat = "identity")+
xlab("\nTreatment\n") + # X axis label
ylab("\nMean Vessicle Concentration\n(10^9 particles/ ml)\n") + # Y axis label
ggtitle("Effect of BPS on Extracellular Vessicle\nRelease of hCTBs\n")+
scale_y_continuous(breaks = seq(0,14,2),
limits = c(0,14),
expand = c(0,0))+ # set bottom of graph
labs(fill= "Time (hr)")
plot2
# ggsave("BPS_treated_hCTBs_48_96_plot.png",
# height = 5, width = 7, dpi = 600, units= "in")Interactive Plot
ggplotly(plot2)Statistics
fit <- aov(corrected_particle_conc ~ Time * Treatment ,data=merge3)
tidy(fit)## term df sumsq meansq
## 1 Time 1 13.57660 13.576598
## 2 Treatment 3 21.35316 7.117719
## 3 Time:Treatment 3 32.29186 10.763954